Understanding patterns of species abundance and persistence in TROLL
In this document, I examine how a species’ performance in mixture depends on 1) its monoculture biomass and 2) its trait values. To measure performance in mixture, I use two indices: 1) total biomass and 2) time to extinction.
1 How trait values translate into monoculture biomass
First, a facetted plot of the relationship between each species’ trait values and its corresponding monoculture biomass.
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year %in% c(100, 200)) %>%
rename(monoculture.biomass = Biomass) %>%
select(-Ninitial, -Year, -Rep) %>%
pivot_longer(c(-SpeciesID, -Stage, -monoculture.biomass),
names_to = "trait",
values_to = "trait.value")
troll.monocultures <- troll.monocultures %>%
group_by(trait) %>%
nest() %>%
summarise(plot = map(.x = data,
.f = ~ {
ggplot(.x) +
geom_point(aes(x = trait.value,
y = monoculture.biomass,
color = SpeciesID)) +
facet_grid(cols = vars(Stage)) +
labs(x = "Trait value",
y = "Monoculture biomass") +
theme_bw() +
theme(aspect.ratio = 1,
legend.position = "none")
}))
cowplot::plot_grid(plotlist = troll.monocultures$plot,
labels = troll.monocultures$trait)
Color indicates SpeciesID. From visual inspection, it appears that dmax is the most important trait determining monoculture biomass.
Because I’m going to be plotting the species’ traits in relation to each other and their species’ emergent biomass, I need to pick which traits are the most important. Therefore, I use a random forest to find the importance for each of the traits in predicting monoculture biomass. Note that, to increase the amount of data for the random forest, I use the last 10 years of the simulation as input data for the random forest, adding Year as a predictor to ensure that temporal autocorrelation does not bias the variable importance measures (see here).
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year %in% seq(180, 200)) %>%
rename(monoculture.biomass = Biomass) %>%
mutate(id = row_number()) %>%
mutate_if(is.character, as.factor) %>%
select(-SpeciesID, -Ninitial, -Stage, -Rep)
train <- troll.monocultures %>% sample_frac(.70)
test <- anti_join(troll.monocultures, train, by = 'id')
train <- train %>% select(-id)
test <- test %>% select(-id)
rf <- randomForest(data = as.data.frame(train),
monoculture.biomass ~ .,
importance = TRUE)
pred <- data.frame(pred = predict(rf, test))
title = paste("correlation: ", round(cor(pred, test$monoculture.biomass)[[1]], 2),
" | mean squared error: ", round(mean(rf$mse), 2),
" | R-squared: ", round(mean(rf$rsq), 2),
sep = "")
varImpPlot(rf, main = title)
Given these importance values, I use dmax and wsg as the focal traits.
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year == 200) %>%
select(-Ninitial, -Stage, -Year, -Rep)
ggplot(troll.monocultures) +
geom_point(aes(x = dmax,
y = wsg,
color = Biomass),
size = 4) +
scale_color_viridis() +
labs(color = "Monoculture\nbiomass") +
theme_bw()
High monoculture biomass is generally related to having high wsg and dmax values.
2 Competitive ability in mixture vs. monoculture biomass
I now explore how the competitiveness of a species in a community (either 32- or, later, 64-species) depends on 1) its monoculture biomass, and subsequently 2) its trait values.
2.1 32-species mixtures
2.1.1 Species’ biomass in 32-species mixture vs. its biomass in monoculture
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year == 200) %>%
rename(monoculture.biomass = Biomass) %>%
select(-Ninitial, -Stage, -Year, -Rep)
troll.32 <- troll %>%
filter(Ninitial == 32,
Year == 200) %>%
select(SpeciesID, Biomass, dmax, wsg, lma, pmass, nmass, h_realmax)
troll.combined <- inner_join(troll.monocultures, troll.32)
ggplot(troll.combined) +
geom_point(aes(x = monoculture.biomass,
y = Biomass,
color = SpeciesID)) +
labs(x = "Monoculture biomass",
y = "Biomass in 32-species mixture") +
theme_bw() +
theme(legend.position = "none")
Each color represents a species, while the number of discrete points per species depends on the number of replicates the species was included in. Inspecting this relationship visually, it appears that species with low biomass monocultures stand a better chance to dominate 32-species mixtures. However, many relatively high biomass monocultures remain in the community.
2.1.2 Biomass in 32-species mixture vs. trait values
To understand this relationship more deeply, I inspect how species’ traits impact their performance in 32-species mixtures.
troll.32 <- troll %>%
filter(Ninitial == 32,
Year == 200) %>%
select(Stage, SpeciesID, Biomass, dmax, wsg, lma, pmass, nmass, h_realmax)
troll.32 <- troll.32 %>%
pivot_longer(c(-Stage, -SpeciesID, -Biomass),
names_to = "trait",
values_to = "trait.value")
troll.32 <- troll.32 %>%
group_by(trait) %>%
nest() %>%
summarise(plot = map(.x = data,
.f = ~ {
ggplot(.x) +
geom_point(aes(x = trait.value,
y = Biomass,
color = SpeciesID)) +
facet_grid(cols = vars(Stage)) +
labs(x = "Trait value",
y = "Monoculture biomass") +
theme_bw() +
theme(aspect.ratio = 1,
legend.position = "none")
}))
cowplot::plot_grid(plotlist = troll.32$plot,
labels = troll.32$trait,
ncol = 2)
It is difficult to discern which trait is the most important predictor of total biomass within the 32-species mixtures. Therefore, I again employ a randomForest to find the traits’ importances.
To increase the data available to the random forest, I again use the last 20 years of the dataset to generate the random forest model, incorporating Year as a predictor to ensure that temporal autocorrelation plays a minimal role in biasing the results.
troll.32 <- troll %>%
filter(Ninitial == 32,
Year %in% seq(180, 200))
troll.32 <- troll.32 %>%
ungroup() %>%
mutate(id = row_number()) %>%
mutate_if(is.character, as.factor) %>%
select(-Ninitial, -Rep, -Stage, -SpeciesID)
train <- troll.32 %>% sample_frac(.70)
test <- anti_join(troll.32, train, by = 'id')
train <- train %>% select(-id)
test <- test %>% select(-id)
rf <- randomForest(data = as.data.frame(train),
Biomass ~ .,
importance = TRUE)
pred <- data.frame(pred = predict(rf, test))
title = paste("correlation: ", round(cor(pred, test$Biomass)[[1]], 2),
" | mean squared error: ", round(mean(rf$mse), 2),
" | R-squared: ", round(mean(rf$rsq), 2),
sep = "")
varImpPlot(rf, main = title)
lma is the best predictor of total biomass in mixture, followed by dmax and pmass, though all traits play an important roll in determining a species’ biomass within a community. Year is a poor predictor of species’ biomass in monoculture, indicating that sampling the last 20 years doesn’t skew the importance values.
A striking result is:
The traits predicting monoculture biomass are not consistent with those that predict biomass in mixture. This corroborates the observation that monoculture biomass does not correlate cleanly with biomass in mixture.
troll.32 <- troll %>%
filter(Ninitial == 32,
Year == 200) %>%
group_by(SpeciesID, dmax, wsg, lma, pmass, nmass, h_realmax) %>%
summarise(Biomass = mean(Biomass)) %>%
filter(Biomass > 0)
plot_ly(troll.32,
x = ~dmax,
y = ~pmass,
z = ~lma,
color = ~log(Biomass),
mode = "markers",
type = "scatter3d")Generally speaking, species with low dmax, low pmass, and low lma dominate the mixtures at the end of the experimental period. However, species with high dmax values are also present to a lesser degree.
2.1.3 Time to extinction in 32-species mixture vs. biomass in monoculture
I now invesigate another metric of competition, time to extinction.
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year == 200) %>%
rename(monoculture.biomass = Biomass) %>%
select(-Ninitial, -Stage, -Year, -Rep)
troll.32 <- troll %>%
filter(Ninitial == 32)
troll.32 <- troll.32 %>%
group_by(SpeciesID, Rep) %>%
arrange(Year) %>%
mutate(isDead = (Biomass == 0)) %>%
filter(isDead | (!isDead & Year == max(Year))) %>%
filter(Year == first(Year))
troll.32 <- troll.32 %>%
group_by(SpeciesID) %>%
summarise(Year = mean(Year))
troll.combined <- inner_join(troll.monocultures, troll.32)
ggplot(troll.combined) +
geom_point(aes(x = monoculture.biomass,
y = Year,
color = SpeciesID)) +
labs(x = "Monoculture biomass",
y = "Extinction year") +
theme_bw() +
theme(legend.position = "none")
Generally speaking, the species surviving to the end of the simulation tend to have lower monoculture biomasses. However, this is not insured; many low monoculture biomass species go extinct quickly.
2.1.4 Time to extinction in 32-species mixture vs. traits
troll.32 <- troll %>%
filter(Ninitial == 32)
troll.32 <- troll.32 %>%
group_by(SpeciesID, Rep) %>%
arrange(Year) %>%
mutate(isDead = (Biomass == 0)) %>%
filter(isDead | (!isDead & Year == max(Year))) %>%
filter(Year == first(Year))
troll.32 <- troll.32 %>%
group_by(SpeciesID, dmax, wsg, lma, pmass, nmass, h_realmax) %>%
summarise(Year = mean(Year))
troll.32 <- troll.32 %>%
pivot_longer(c(-SpeciesID, -Year),
names_to = "trait",
values_to = "trait.value")
ggplot(troll.32) +
geom_point(aes(x = trait.value,
y = Year,
color = SpeciesID)) +
facet_wrap(facets = vars(trait),
ncol = 4,
scales = "free") +
labs(x = "Trait value",
y = "Time to extinction") +
theme_bw() +
theme(legend.position = "none")
From this graph, we can see that the trait lma is likely the strongest factor predicting time to extinction. To check this intuition, I again use a random forest.
troll.32 <- troll %>%
filter(Ninitial == 32)
troll.32 <- troll.32 %>%
group_by(SpeciesID, Rep) %>%
arrange(Year) %>%
mutate(isDead = (Biomass == 0)) %>%
filter(isDead | (!isDead & Year == max(Year))) %>%
filter(Year == first(Year)) %>%
rename(extinction.year = Year)
troll.32 <- troll.32 %>%
ungroup() %>%
mutate(id = row_number()) %>%
mutate_if(is.character, as.factor) %>%
select(-SpeciesID, -Ninitial, -Stage, -Rep, -Biomass, -isDead)
train <- troll.32 %>% sample_frac(.70)
test <- anti_join(troll.32, train, by = 'id')
train <- train %>% select(-id)
test <- test %>% select(-id)
rf <- randomForest(data = as.data.frame(train),
extinction.year ~ .,
importance = TRUE)
pred <- data.frame(pred = predict(rf, test))
title = paste("correlation: ", round(cor(pred, test$extinction.year)[[1]], 2),
" | mean squared error: ", round(mean(rf$mse), 2),
" | R-squared: ", round(mean(rf$rsq), 2),
sep = "")
varImpPlot(rf, main = title)
Confirming my intuition, lma is by far the best predictor of time to extinction, followed by pmass.
troll.32 <- troll %>%
filter(Ninitial == 32)
troll.32 <- troll.32 %>%
group_by(SpeciesID, Rep) %>%
arrange(Year) %>%
mutate(isDead = (Biomass == 0)) %>%
filter(isDead | (!isDead & Year == max(Year))) %>%
filter(Year == first(Year))
troll.32 <- troll.32 %>%
group_by(SpeciesID, dmax, wsg, lma, pmass, nmass, h_realmax) %>%
summarise(Year = mean(Year))
ggplot(troll.32) +
geom_point(aes(x = lma,
y = pmass,
color = Year)) +
labs(color = "Extinction\nyear") +
scale_color_viridis() +
theme_bw()
Species with low lma and pmass values remain in the community for a dramatically longer time species with high values.
A potential insight:
When considering how well adapted a species is to persist within these communities, trait values that convey persistence (time to extinction) are less complex to interpet than those converying abundance.
2.2 64-species mixture
2.2.1 Biomass in the 64-species mixture vs. monoculture biomass
To understand how these dynamics change with seed addition, I incorporate the metacommunity phase into the following analysis of 64-species mixtures. During the metacommunity phase, seed addition is present.
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year %in% c(100, 200)) %>%
rename(monoculture.biomass = Biomass) %>%
select(-Ninitial, -Year, -Rep)
troll.mixture <- troll %>%
filter(Ninitial == 64,
Year %in% c(100, 200)) %>%
rename(mixture.biomass = Biomass) %>%
select(-Ninitial, -Year, -Rep)
troll.combined <- inner_join(troll.monocultures, troll.mixture)
ggplot(troll.combined) +
geom_point(aes(x = monoculture.biomass,
y = mixture.biomass,
color = SpeciesID)) +
facet_grid(cols = vars(Stage)) +
labs(x = "Monoculture biomass",
y = "Mixture biomass") +
theme_bw() +
theme(legend.position = "none")
It appears that there is a switch between the two phases of the experiment. External seed addition tends to prop up the high monoculture biomass species such that they dominate the mixtures as well. However, after seed addition is removed, species with low monoculture biomasses tend to dominate.
2.2.2 Time to extinction in the 64-species mixture vs. monoculture biomass
troll.monocultures <- troll %>%
filter(Ninitial == 1,
Year == 200) %>%
rename(monoculture.biomass = Biomass) %>%
select(-Ninitial, -Year, -Rep)
troll.64 <- troll %>%
filter(Ninitial == 64) %>%
select(-Ninitial, -Rep)
troll.64 <- troll.64 %>%
group_by(SpeciesID) %>%
arrange(Year) %>%
mutate(isDead = (Biomass == 0)) %>%
filter(isDead | (!isDead & Year == max(Year))) %>%
filter(Year == first(Year))
troll.combined <- inner_join(troll.monocultures, troll.64)
ggplot(troll.combined) +
geom_point(aes(x = monoculture.biomass,
y = Year,
color = SpeciesID)) +
labs(x = "Monoculture biomass",
y = "Extinction year") +
theme_bw() +
theme(legend.position = "none")
Because there is only one 64-species mixture, and thus a limited amount of data, I refrain from running a random forest of the results to understand the trait analysis.
3 Conclusions
While dmax and wsg are the most important traits in determining the monoculture biomass of a species within TROLL, lma is by far more important in predicting both the biomass in mixture of a species, as well as how long it can persist within a community.